/**
 * 
 */
package edu.usc.epigenome.dnaase.script;

import java.io.BufferedReader;
import java.io.File;
import java.io.FileReader;
import java.io.IOException;
import java.io.PrintWriter;
import java.util.ArrayList;
import java.util.HashMap;
import java.util.List;

import net.sf.samtools.SAMFileReader;
import net.sf.samtools.SAMRecord;
import net.sf.samtools.util.CloseableIterator;

import org.kohsuke.args4j.Argument;
import org.kohsuke.args4j.CmdLineException;
import org.kohsuke.args4j.CmdLineParser;
import org.kohsuke.args4j.Option;

/**
 * @author yaping
 * @contact lyping1986@gmail.com
 * @time May 6, 2012 9:01:42 PM
 * 
 */
public class GetVCFposHaveCoverage {
	final private static String USAGE = "GetVCFposHaveCoverage [opts] bamFileName oldVcfFileName newVcfFileName";
	/**
	 * @param args
	 */
	@Option(name="-minCov",usage="minimum covergae required for the position in VCF file, default:1")
	protected int minCov = 1;
	@Option(name="-maxCov",usage="maximum covergae required for the position in VCF file, default:120")
	protected int maxCov = 120;
	@Option(name="-minMapQ",usage="Minimum mapping quality in the covered loci when counting coverage, default:30")
	protected int minMapQ = 30;
	@Option(name="-minBaseQ",usage="Minimum base quality in the covered loci when counting coverage, default:0")
	protected int minBaseQ = 0;
	
	@Argument
	private List<String> stringArgs = new ArrayList<String>();
	
	public static void main(String[] args) throws Exception {
		// TODO Auto-generated method stub
		new GetVCFposHaveCoverage().doMain(args);
	}

	public void doMain(String[] args) throws IOException{
		CmdLineParser parser = new CmdLineParser(this);
		// if you have a wider console, you could increase the value;
		// here 80 is also the default
		parser.setUsageWidth(80);
		try
		{
			parser.parseArgument(args);
			if (stringArgs.size() < 3) throw new CmdLineException(USAGE);
		
		}
		catch (CmdLineException e)
		{
			System.err.println(e.getMessage());
			// print the list of available options
			parser.printUsage(System.err);
			System.err.println();
			return;
		}

		
		String bamFileName = stringArgs.remove(0);
		String vcfFileNameInput = stringArgs.remove(0);
		String vcfFileNameOutput = stringArgs.remove(0);
		BufferedReader br = new BufferedReader(new FileReader(vcfFileNameInput));
		SAMFileReader inputSam = new SAMFileReader(new File(bamFileName));
		inputSam.setValidationStringency(SAMFileReader.ValidationStringency.SILENT);
		
		PrintWriter outVcfWriter = new PrintWriter(new File(vcfFileNameOutput));
		int recCounter = 1;
		int oldRecord = 0;
		int newRecord = 0;
		String line;
		while( (line = br.readLine()) != null){
			if(line.startsWith("#")){
				outVcfWriter.println(line);
			}
			else{
				oldRecord++;
				String[] tmpArray = line.split("\t");
				String chr = tmpArray[0];
				int location = Integer.parseInt(tmpArray[1]);
				CloseableIterator<SAMRecord> chrIt = inputSam.queryOverlapping(chr, location, location);
				int converage = 0;
				while (chrIt.hasNext())
				{	
					SAMRecord samRecord = chrIt.next();
					if(samRecord.getDuplicateReadFlag() || samRecord.getNotPrimaryAlignmentFlag() || samRecord.getReadFailsVendorQualityCheckFlag() || samRecord.getReadUnmappedFlag() || samRecord.getMappingQuality() <= minMapQ)
						continue;
					if(samRecord.getReadPairedFlag()){
						if (samRecord.getAlignmentStart() == samRecord.getMateAlignmentStart() && samRecord.getReadNegativeStrandFlag() == samRecord.getMateNegativeStrandFlag()){
							if (samRecord.getSecondOfPairFlag()) 
		 						continue;
						}
		 					if(!samRecord.getProperPairFlag())
		 						continue;
					}
					converage++;	
				}
				chrIt.close();
				
				if(converage>=minCov && converage<=maxCov){
					outVcfWriter.println(line);
					newRecord++;
				}
				if ((recCounter % 30000)==0)
				{
					System.err.printf("On new record #%d\n",recCounter);
					inputSam.close();
					
					System.gc();
					inputSam = new SAMFileReader(new File(bamFileName));
					inputSam.setValidationStringency(SAMFileReader.ValidationStringency.SILENT);
					
				}
			}
			
		}
		br.close();
		inputSam.close();
		outVcfWriter.close();
		System.out.println("finished!");
		System.out.println("old VCF record: " + oldRecord);
		System.out.println("old VCF record: " + newRecord);
	}
}
